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' The genetic cluster-exact approximation algorithm is an 

efficient method to calculate ground states of EA spin glasses. 
The method can be used to study ground-state landscapes by 
calculating many independent ground states for each realiza- 
tion of the disorder. The algorithm is analyzed with respect 
to the statistics of the ground states and the valleys of the 
energy landscape. Furthermore, the distribution inside each 
valley is evaluated. It is shown that the algorithm does not 
lead to a true T — thermodynamic distribution, i.e. each 
ground state has not the same frequency of occurrence when 
performing many runs. An extension of the technique is out- 
lined, which guarantees that each ground states occurs with 
the same probability. 

Keywords (PACS-codes): Spin glasses and other 
random models (75.10.Nr), Numerical simulation studies 
(75.40.Mg), General mathematical systems (02.10.Jf). 



I. INTRODUCTION 

The finite-dimensional Edwards-Anderson spin glass 
is a model for disordered systems which has attracted 
much attention over the last decades. The opinion on its 
nature, especially for three dimensional systems, is still 
controversial [^-Q . Beside trying to address the problem 
with the help of analytic calculations and simulations at 
finite temperature, it is possible to investigate the behav- 
ior of the model by means of ground-state calculations . 
Since obtaining spin-glass ground states is computation- 
ally hard the study is restricted to relatively small 
systems. Recently a new algorithm, the cluster-exact ap- 
proximation (CEA) was presented, which allows in 
connection with a special genetic algorithm the cal- 
culation of true ground states for moderate system 
sizes, in three dimensions up to size 14'^. By applying this 
method it is possible to study the ground-state landscape 
of systems exhibiting a T = degeneracy . For a ther- 
modynamical correct evaluation it is necessary that each 
ground state contributes to the results with the same 
weight, since all ground states have exactly the same en- 
ergy. Recently it was shown [Q, that the genetic CEA 
causes a bias on the quantities describing the T = land- 
scape. The aim of this paper is to analyze the algorithm 
with respect to its ground-state statistics. The reasons 
for the deviation from the correct behavior are given and 



an extension of the method is outlined, which guarantees 
thermodynamical correct results. 

In this work, three-dimensional Edwards-Anderson 
(EA) ±J spin glasses are investigated. They consist of 
N spins CTi — ±1, described by the Hamiltonian 
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The sum runs over all pairs of nearest neighbors. The 
spins are placed on a three-dimensional (d=3) cubic lat- 
tice of linear size L with periodic boundary conditions in 
all directions. Systems with quenched disorder of the in- 
teractions (bonds) are considered. Their possible values 
are = ±1 with equal probability. To reduce the fluc- 
tuations, a constraint is imposed, so that j) ^ij — ^■ 

The article is organized as follows: next a description 
of the algorithms is presented. Then it is shown for small 
systems, that the method does not result in a thermo- 
dynamical correct distribution of the ground states. In 
section four, the algorithm and its different variants are 
analyzed with respect to the ground-state statistics. In 
the last section a summary is given and an extension of 
the method is outlined, which should guarantee thermo- 
dynamical correct results. 



II. ALGORITHMS 

The algorithm for the calculation bases on a special 
genetic algorithm ||l^,|l^ and on cluster-exact approxi- 
mation |10|. CEA is an optimization method designed 
specially for spin glasses. Its basic idea is to transform 
the spin glass in a way that graph-theoretical methods 
can be applied, which work only for systems exhibiting 
no bond-frustrations. Now a short sketch of these algo- 
rithms is given, because later the influence of different 
variants on the results is discussed. 

Genetic algorithms are biologically motivated. An op- 
timal solution is found by treating many instances of the 
problem in parallel, keeping only better instances and 
replacing bad ones by new ones (survival of the fittest). 
The genetic algorithm starts with an initial population 
of Mi randomly initialized spin configurations (= indi- 
viduals) , which are linearly arranged using an array. The 
last one is also neighbor of the first one. Then rio x Mi 
times two neighbors from the population are taken (called 
parents) and two new configurations called offspring are 
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created. For that purpose the triadic crossover is used 
which turned out to be very efficient for spin glasses: a 
mask is used which is a third randomly chosen (usually 
distant) member of the population with a fraction of 0.1 
of its spins reversed. In a first step the offspring are 
created as copies of the parents. Then those spins are se- 
lected, where the orientations of the first parent and the 
mask agree |16 . The values of these spins are swapped 



between the two offspring. Then a mutation with a rate 
of Pm is applied to each offspring, i.e. a randomly chosen 
fraction of the spins is reversed. 

Next for both offspring the energy is reduced by ap- 
plying CEA: The method constructs iteratively and ran- 
domly a non- frustrated cluster of spins. During the con- 
struction of the cluster a local gauge-transformation of 
the spin variables is applied so that all interactions be- 
tween cluster spins become ferromagnetic. Fig. |^ shows 
an example of how the construction of the cluster works 
for a small spin-glass system. To increase the perfor- 
mance, spins adjacent to many unsatisfied bonds are 
more likely to be added to the cluster. This may intro- 
duce a bias on the resulting distribution of the ground 
states. Later this scheme ("BIAS") is compared to a 
variant ("SAME"), where all spins may contribute to the 
cluster with the same probability. 

For 3d ± J spin glasses each cluster contains typically 
55 percent of all spins. The non-cluster spins remain 
fixed during the following calculation, they act like local 
magnetic fields on the cluster spins. Consequently, the 
ground state of the gauge-transformed cluster is not triv- 
ial, although all interactions inside the cluster are ferro- 
magnetic. Since the cluster exhibits no bond-frustration, 
an energetic minimum state for its spins can be calculated 
in polynomial time by using graph-theoretical methods 
p7|-p^: an equivalent network is constructed the 
maximum fiow is calculated ||2l|,^ and the spins of the 
cluster are set to orientations leading to a minimum in 
energy. Please note, that the ground state of the clus- 
ter is degenerate itself, i.e. the spin orientations can be 
chosen in different ways leading all to the same energy. 
It is possible to calculate within one single run a special 
graph, which represents all ground states of the cluster 
|p3| , and select one ground state randomly. This proce- 
dure is called "BROAD" here. On the other hand, one 
can always choose a certain ground state of the cluster 
directly]]. Usually this variant, which is called "QUICK" 
here, is applied, because it avoids the construction of the 
special graph. But this again introduces a certain bias 
on the resulting distribution of the ground states. Later 
the influence of the different methods of choosing ground 
states is discussed. 



^This ground state has the maximum possible magnetiza- 
tion of the gauge-transformed spins among all cluster ground 
states. 



This CEA minimization step is performed n„iin times 
for each offspring. Afterwards each offspring is compared 
with one of its parents. The offspring/parent pairs are 
chosen in the way that the sum of the phenotypic differ- 
ences between them is minimal. The phenotypic differ- 
ence is defined here as the number of spins where the two 
configurations differ. Each parent is replaced if its energy 
is not lower (i.e. not better) than the corresponding off- 
spring. After this whole step is conducted rio x Mi times, 
the population is halved: From each pair of neighbors 
the configuration which has the higher energy is elimi- 
nated. If more than 4 individuals remain the process is 
continued otherwise it is stopped and the best individual 
is taken as result of the calculation. 

The following representation summarizes the algo- 
rithm. 

algorithm genetic CEA({ Jij}, Mi, Uo, Pm, n-min) 
begin 

create Mi configurations randomly 

while {M^ > 4) do 

begin 

for i = 1 to Ho X Mi do 

begin 

select two neighbors 

create two offspring using triadic crossover 
do mutations with rate pm 
for both offspring do 
begin 

for j = 1 to n,nin do 

begin 

construct unfrustrated cluster of spins 
construct equivalent network 
calculate maximum flow 
construct minimum cut 
set new orientations of cluster spins 
end 

if offspring is not worse than related parent 
then 

replace parent with offspring 

end 
end 

half population; M^ = Mi/2 
end 

return one configuration with lowest energy 
end 

The whole algorithm is performed rifj times and all 
configurations which exhibit the lowest energy are stored, 
resulting in no statistically independent ground-state 
configurations (replicas). A priori nothing about the 
distribution of ground states raised by the algorithm is 
known. Thus, it may be possible that for one given real- 
ization of the disorder some ground states are more likely 
to be returned by the procedure than others. Conse- 
quently, any quantities which are calculated by averaging 
over many independent ground states, like the distribu- 
tion of overlaps, may depend on a bias introduced by the 
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algorithm. For a thermodynamical correct evaluation all 
ground states have to contribute with the same weight, 
since they all have exactly the same energy. 

For the preceding work, the distribution of the ground 
states determined by the algorithm was taken. The 
method was utilized to examine the ground state land- 
scape of two-dimensional p4| ] and three-dimensional 
|^ , ^5| ±J spin glasses by calculating a small number 
of ground states per realization. Some of these results 
depend on the statistics of the ground states, as it will 
be shown in the next section for the d = 3 case. 

On the other hand, the main findings of the follow- 
ing investigations are not affected by the bias intro- 
duced by genetic CEA: the existence of a spin-glass phase 
for nonzero temperature was confirmed for the three- 
dimensional spin glass [ p^ . The method was applied also 
to the ± J random-bond model to investigate its T = 
ferromagnetic to spin-glass transition [ p6| . Finally, for 
small sizes up to L = 8 all ground-state valleys were 
obtained by calculating a huge number of ground states 
per realization and applying a new method called ballistic 
search 

III. NUMERICAL EVIDENCE 

In this section results describing the ground-state land- 
scape of small three-dimensional ± J spin glasses are eval- 
uated. It is shown that the data emerging from the use of 
raw genetic CEA and from a thermodynamically correct 
treatment differ substantially. 

Several ground states for small systems of size iV = 
= 3"^, 4'^, 5'^ were calculated. 1000 realizations of the 
disorder for L = 3, 4 and 100 realizations for L = 5 
were considered. The parameters (M^, rto, p™, Jimin), 
for which true ground states are obtained, are shown in 
p2[ . For all calculations the variants BIAS and QUICK 
were used to obtain maximum performance. The effect 
of different variants on the results is discussed in the next 
section. 

Two schemes of calculation were applied: 

A For each realization nj^ — 40 runs of genetic 
CEA were performed and all states exhibiting the 
ground-state energy stored. Consequently, this 
scheme reflects the ground-state statistics which 
is determined solely by the genetic CEA method. 
Configurations which have a higher probability of 
occurrence contribute with a larger weight to the 
results. 

B For each realization the algorithm was run up to 
10^ times. Each particular state was stored only 
once. For later analysis the number of times each 
state occurred was recorded. Additionally, a sys- 
tematic local search was applied to add possibly 
missing ground states which are related by flips of 
free spins to states already found. Finally, a L = 3 



realization exhibits 25 different ground states on 
average. For a L = 4 realization on average 240 
states were found and 6900 states for L — 5. 

For the evaluation of physical quantities every 
ground state is taken with the same probability in 
this scheme. Thus, the statistics obtained in this 
way reflect the true T — thermodynamic behav- 
ior. 

To analyze the ground-state landscape, the distribu- 
tion of overlaps is evaluated. For a fixed realization 
J — {Jij} of the exchange interactions and two replicas 
{(y'i^ , {<y^^ , the overlap ||^ is defined as 

'^"'-^E-r-f (2) 

i 

The ground state of a given realization is characterized 
by the probability density Pj{q). Averaging over the 
realizations J, denoted by [•],/, results in {Z = number 
of realizations) 

P(g)^[P;(<z)], = l^P,(g) (3) 
J 

Because no external field is present the densities are 
symmetric: Pj{q) = Pj{~q) and P{q) ~ P{—q). So only 
P(|g|) is relevant. 

The result of P(|g|) for L = 5 is shown in Fig. For 
the true thermodynamic result small overlaps occur less 
frequent than for the data obtained by the application 
of pure genetic CEA. Large overlap values occur more 
often. This deviation has an infiuence on the way the 
spin glass behavior is interpreted. The main controversy 
about finite-dimensional spin glasses mentioned at the 
beginning is about the question whether for the infinite 
system P(|<7|) shows a long tail down to g = or not 

To investigate the finite size behavior of P(|g|) the frac- 
tion ATq.s of the distribution below qo = 0.5 is integrated: 

rqo 

X,, ^ / P{\q\)dq (4) 
Jo 

The development of ATq.s as a function of system size 
L is shown in Fig. ^. The datapoints for the larger sizes 
L > 6, obtained using pure genetic CEA, are taken from 
former calculations |]l3|. These values are more or less 
independent of the system size, while the correct thermo- 
dynamic behavior shows a systematic decrease. Whether 
for L ^ oo the long tail of P(|(7|) persists cannot be con- 
cluded from the data, because the systems are too small. 
Nevertheless, the true T = behavior differs significantly 
from the former results. 



IV. ANALYSIS OF GENETIC CEA 

To understand, why genetic CEA fails in producing 
the thermodynamical correct results, in this section the 



3 



statistics of the ground states, which is determined by 
the algorithm, is analyzed directly. 

For the case where all ground states were calculated us- 
ing a huge number of runs, the frequencies each ground 
state occurred were recorded. In Fig. ^ the result for 
one sample realization of = 5"^ is shown. The system 
has 56 different ground states. For each state the num- 
ber of times it was returned by the algorithm in 10^ runs 
is displayed. Obviously the large deviations from state 
to state cannot be explained by the presence of statis- 
tical fluctuations. Thus, genetic CEA samples different 
ground states from the same realization with different 
weights. 

To make this statement more precise, the following 
analysis was performed: Two ground states are called 
neighbors, if they differ only by the orientation of one 
spin. All ground states which are accessible from each 
other through this neighbor-relation are defined to be in 
the same ground-state valley. That means, two ground 
states belong to the same valley, if it is possible to move 
from one state to the other by flipping only free spins, 
i.e. without changing the energy. For all realizations 
the valleys were determined using a method presented in 
, which allows to treat systems efficiently exhibiting a 
huge number of ground states. Then the frequencies hv 
for each valley V were computed as the sum of all fre- 
quencies of the states belonging to V. In Fig. ||the result 
is shown for a sample = 5'^ realization, which has 15 
different ground state valleys. Large valleys are returned 
by the algorithm more frequently, but hy seems to grow 
slower than linearly. A strict linear behavior should hold 
for an algorithm which guarantees the correct T = be- 
havior. 

For averaging hv has to be normalized, because the 
absolute values of the frequency differ strongly from re- 
alization to realization, even if the size \V\ of a valley, 
i.e. the number of ground states belonging to it, is the 
same. For each realization, the normalized frequency hy 
is measured relatively to the average frequency hi of all 
valleys of size 1: hy = hy /hi 

If a realization does not exhibit a valley consisting only 
of one ground state, the frequency hy^ of the smallest 
valley Vg is taken. It is assumed, that the normalized 
frequency exhibits a hy — dependence, which is 

justified by the results shown later. Consequently, for the 
case the size \Vs \ of the smallest valley is larger than one, 
hi = hv^/iVsl" is chosen. The value of a is determined 
self-consistently. 

The result for L — 3 of /ly as a function of the valley- 
size \V\ is presented in Fig. |6[ A value of a = 0.854(3) 
was determined. Please note, that the fluctuations for 
larger valleys are higher, because quite often only one 
valley was available for a given valley-size. The algebraic 
form is clearly visible, proving that genetic CEA overes- 
timates systematically the importance of small ground- 
state valleys. 

For L = 4 a value oi a — 0.705(3) was obtained. 



while the L = 5 case resulted in a = 0.642(5). Conse- 
quently, with increasing system size, the algorithm fails 
more and more to sample configurations from different 
ground-state valleys according to the size of the valleys. 
This explains, why the difference of Xo,5{L) between the 
correct result and the values obtained in increases 
with growing system size. 

Similar results were obtained for two-dimensional sys- 
tems. For i = 5 a self-consistent value of a = 0.650(1) 
was found, while the treatment of L = 7 systems resulted 
in a = 0.659(2). Here only a slight finite-size dependence 
occurs. This may explain the fact, that the width of the 
distribution of overlaps, even calculated only by the ap- 
plication of pure genetic CEA, seems to scale to zero [ p4[ . 

In the second section of this paper two variants of the 
algorithm were presented, which may be able to calcu- 
late ground states more equally distributed. To inves- 
tigate this issue, similar ground-state calculations were 
conducted for L — 4 and again hy was calculated. For 
the case, were SAME was used instead of BIAS, a value 
a = 0.801(2) was determined self-consistently. Using 
BROAD instead of QUICK resulted in a = 0.749(3). 
Finally, by applying SAME and BROAD together, a = 
0.843(3) was obtained. Consequently, applying different 
variants of the method decreases the tendency of over- 
estimating small valleys, but the correct thermodynamic 
behavior is not obtained as well. Even worse, BROAD 
and SAME are considerably slower than the combination 
of QUICK and BIAS. 

So far it was shown, that genetic CEA fails in sampling 
ground states from different valleys according the size of 
the valleys. Now we turn to the question, whether at least 
states belonging to the same valley are calculated with 
the correct thermodynamic distribution. By investigat- 
ing the frequencies of different ground states belonging 
to the same valley it was found again, that these config- 
urations are not equally distributed. But it is possible 
to study this issue in a more physical way. For that 
purpose ground states of 100 L = 10 realizations were 
calculated. Then the valley structure was analyzed. The 
average distribution of overlaps was evaluated, but only 
contributions of pairs of states belonging to the same 
valley were considered. For comparison, for the same re- 
alizations a long T = Monte-Carlo (MC) simulation 
was performed, i.e. randomly spins were selected and 
flipped if they were free. The ground states were used as 
starting configurations. Since a MC simulation ensures 
the correct thermodynamic distribution of the states, all 
ground states of a valley appear with the same frequency, 
if the simulation is only long enough. A length of 40 
Monte-Carlo steps per spin were found to be sufficient 
for L = 10. The result for the distribution of overlaps 
^valley(l^l) restricted to the valleys is displayed in Fig. 
0. Significant differences between the datapoints from 
the pure genetic CEA and the correct T = behavior 
are visible. Consequently, the algorithm does not sample 
configurations belonging to the same ground-state valley 
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with the same weight as well. 

V. CONCLUSION 

In this work the genetic cluster-exact approximation 
method is analyzed. The algorithm can be used to calcu- 
late many independent true ground states of EA spin 
glasses. The results from the raw application of the 
method and from calculations of all ground states for 
small system sizes were compared. By evaluating the 
distribution of overlaps is was shown, that genetic CEA 
imposes some bias on the ground-state statistics. Con- 
sequently, the results from the application of the raw 
method do not represent the true T = thermodynam- 
ics. 

To elucidate the behavior of the algorithm the statis- 
tics of the ground states were evaluated directly. It was 
shown, that different ground states have dissimilar prob- 
abilities of occurrence. To understand this effect better, 
the ground-state valleys were determined. The genetic 
CEA method finds configurations from small ground- 
state valleys relative to the size of the valley more of- 
ten than configurations from large valleys. Additionally, 
within a valley the states are not sampled with the same 
weight as well. It was shown that two variants of the al- 
gorithm, which decrease its efficiency, weaken the effect, 
but it still persists. 

Summarizing, two effects are responsible for the biased 
ground-state sampling of genetic CEA: small valleys are 
sampled too frequently and the distribution within the 
valleys is not flat. 

For small system sizes it is possible to calculate all 
ground states, so one can obtain the true thermodynamic 
average directly. But already for L = 5 there are realiza- 
tions exhibiting more than 10^ different ground states. 
Since the ground-state degeneracy grows exponentially 
with system size [ p7| larger systems cannot be treated in 
this way. The following receipt should overcome these 
problems and should allow to obtain the true thermody- 
namic T — behavior for larger systems: 

• Calculate several ground states of a realization us- 
ing genetic CEA. 

• Identify the ground states which belong to the same 
valleys. 

• Estimate the size of each valley. This can be done 
using a variant of ballistic search |^ , which works 
by flipping free spins sequentially, each spin at most 
once. The number of spins flipped is a quite accu- 
rate measure for the size of a valley. 

• Sample from each valley a number of ground states, 
which is proportional to the size of the valley. This 
guarantees, that each valley contributes with its 
proper weight. Each state is obtained by per- 
forming a T = MC simulation of sufficient 



length, starting with true ground-state configura- 
tions. Since MC simulations achieve a thermo- 
dynamical correct distribution, it is guaranteed 
that the states within each valley are equally dis- 
tributed. 

Please note, that it is not necessary to calculate all 
ground states to obtain the true thermodynamic behav- 
ior, because it is possible to estimate the size of a valley 
by analyzing only some sample ground states belonging 
to it. Furthermore, it is even only necessary to have con- 
figurations from the largest valleys available, since they 
dominate the ground-state behavior. This condition is 
fulfilled by genetic CEA, because large valleys are sam- 
pled more often than small valleys, even if small valleys 
appear too often relatively. 

From the results presented here it is not possible to 
deduce the correct T = behavior of the infinite system, 
because the system sizes are too small. Using the scheme 
outlined above, it is possible to treat system sizes up to 

L=i4 m. 
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FIG. 1. Example of the cluster-exact approximation 
method. A part of a spin glass is shown. The circles rep- 
resent lattice sites/spins. Straight lines represent ferromag- 
netic bonds the jagged lines antiferromagnetic interactions. 
The top part shows the initial situation. The construction 
starts with the spin at the center. The bottom part displays 
the final stage. The spins which belong to the cluster carry 
a plus or minus sign which indicates how each spin is trans- 
formed, so that only ferromagnetic interactions remain inside 
the cluster. All other spins cannot be added to the cluster 
because it is not possible to multiply them by ±1 to turn all 
adjacent bonds ferromagnetic. Please note that many other 
combinations of spins can be used to build a cluster without 
frustration. 
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FIG. 2. Average distribution of overlaps -P(|?|) for L — 5. 
The dashed line shows the old result obtained by computing 
about 40 independent ground states per realization using ge- 
netic CEA. The solid line shows the same quantity for the 
case, where all existing ground states were used for the evalu- 
ation, i.e. where the correct T = thermodynamic behavior 
is ensured. 
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FIG. 3. Average fraction X0.5 of the distribution of over- 
laps for \q\ < 0.5 as a function of system size L. The upper 
points (circles) where obtained by calculating about 40 inde- 
pendent ground states per realization using genetic CEA. The 
lower points (triangles) show the result (L — 3, 4, 5) for the 
case, where all existing ground states were used for the evalu- 
ation, i.e. where the correct T = thermodynamic behavior 
is ensured. 
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FIG. 4. Number of times each ground state is calculated in 
10^ runs of genetic CEA for a sample realization of A'' = 5^. 
The realization exhibits 56 ground states, which have signifi- 
cant different probabilities of being calculated. 
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FIG. 5. Number of times a ground state belonging to a 
specific valley is calculated in 10^ runs of genetic CEA. The 
result is shown as a function of the valley size \V\ and for one 
L = 5 sample realization. 
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FIG. 6. Normalized number of times a specific valley V is 
found by genetic CEA for L = 3. The frequency is normalized 
so that hy = 1 (see text). The probability that a cluster is 
found increases with the size of the cluster, but slower than 
linearly. The line shows a fit ft;^ = \V\'^ with a = 0.854(3) 




FIG. 7. Distribution ^Valley (I*?!) "-"^ overlaps for L = 10 
restricted to pairs of ground states belonging to the same 
valley. The full line shows the result for the case, where the 
statistics of the ground state is determined by the genetic 
CEA algorithm. The data represented by the dashed line was 
obtained using states which are equally distributed within 
each valley, which was guaranteed by performing a T = MC 
simulation. 
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